Statistical Mechanics of Double sinh-Gordon Kinks 
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We study the classical thermodynamics of the double sinh-Gordon (DSHG) theory in 1+1 
dimensions. This model theory has a double well potential V(</>) = (C cosh 20 — n) 2 when 
n > C, thus allowing for the existence of kinks and antikinks. Though it is nonintegrable, 
the DSHG model is remarkably amenable to analysis. Below we obtain exact single kink and 
kink lattice solutions as well as the asymptotic kink-antikink interaction. In the continuum 
limit, finding the classical partition function is equivalent to solving for the ground state of a 
Schrodinger-like equation obtained via the transfer integral method. For the DSHG model, 
this equation turns out to be quasi-exactly solvable. We exploit this property to obtain 
exact energy eigenvalues and wavefunctions for several temperatures both above and below 
the symmetry breaking transition temperature (provided n = 1, 2, ■ • • , 6). The availability 
of exact results provides an excellent testing ground for large scale Langevin simulations. 
The probability distribution function (PDF) calculated from Langevin dynamics is found to 
be in striking agreement with the exact PDF obtained from the ground state wavefunction. 
This validation points to the utility of a PDF-based computation of thermodynamics utiliz- 
ing Langevin methods. In addition to the PDF, field-field and field fluctuation correlation 
functions were computed and also found to be in excellent agreement with the exact results. 
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PACS numbers: 05.20.-y, ll.10.-z, 63.75,+z, 64.60.Cn 



I. INTRODUCTION 

The double- well </> 4 model in 1 + 1 dimensions has been 
extensively studied in the context of symmetry breaking 
transitions. As is well-known, the equilibrium classical 
statistical mechanics of a 1 + 1 dimensional field theory 
reduces to a time-independent quantum mechanics prob- 
lem via the transfer integral method |l|,||. However, the 
Schrodinger equation with a (ft 4 potential does not have 
any known exact solutions. To overcome this problem, we 
turn to another double-well system that is quasi-exactly 
solvable (QES), where by a QES problem we mean one 
in which the exact partial diagonalization of the Hamilto- 
nian may be carried out [|J. Knowing just the low- lying 
eigenvalues and eigenfunctions is sufficient to obtain al- 
most complete information regarding the classical ther- 
modynamics of the field theoretic system. As will be seen 
below the QES property holds only at a discrete set of 
temperatures. Our strategy is to tune system parame- 
ters so as to straddle an interesting region in tempera- 
ture space which, for this paper, we have taken to be the 
region around the short-range order ("kink") transition 
point. 

An important motivation for exact results is validation 
of numerical techniques: recent advances in large scale 
computation have made it possible to calculate quantities 
such as the probability distribution function (PDF) us- 
ing Langevin dynamics to very high accuracy [IB. The 
PDF is the probability distribution of field values aver- 
aged over the total system volume. Alternatively, it is 



equal to the square of the ground state wavefunction of 
the transfer operator at a given temperature. All ther- 
modynamic information can be shown to reside in the 
PDF J6|,f7| . Calibration against the PDF calculated from 
the exact ground state wavefunction is an essential guide 
when performing error and convergence analysis of the 
Langevin equation la]. 

The double-well QES system that we focus on in this 
paper is specified by the double sinh-Gordon potential 
V(4>) — (C cosh 2(f)— n) 2 , where C is a positive parameter 
and n is a positive integer [|5||j|,[l0 11 l^ , |L3| , |l4| . This po- 



tential is the hyperbolic analog of the double sine-Gordon 
system |l5|] . Similar potentials arise in the context of 
the quantum theory of molecules (e.g., a homonuclear 
diatomic molecule), wave motion describing the normal 
modes of vibration of a stretched membrane of variable 
density [[l6| , and as the solution of a Fokker-Planck equa- 
tion Note that the hyperbolic analog of the sine- 
Gordon equation is a single well potential (sinh-Gordon 
Ip^l ) and thus uninteresting from the soliton statistical 
mechanics perspective. 

Below we find exact solutions for the ground and a few 
excited state energy eigenvalues and wavefunctions for 
certain temperatures both above and below the short- 
range transition temperature. This allows analytic cal- 
culation of the PDF and correlation functions (Cj and 
C2) and direct comparison with those calculated from 
a Langevin simulation. The analytic kink profile can 
also be compared with the corresponding field config- 
uration in the simulation. We find striking agreement 
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between the exact results and simulations. The proven 
high quality of the simulations provides an important val- 
idation for PDF-based thermodynamics implemented via 
Langevin methods. 

Before proceeding to a more detailed view of the DSHG 
model we contrast its main features with those of other 
models that have been used previously for analytic stud- 
ies. The Schrodinger equation with a double-quadratic 
potential V{4>) = \{\<f>\ - l) 2 HQ is in fact exactly 
solvable. However, the energy eigenvalues are not known 
as simple functions but may be obtained only as a so- 
lution of transcendental equations involving parabolic 
cylinder functions. Moreover, this potential is not a 
smooth function and has a cusp at = 0. Other 
smoothly varying double- well potentials, e.g., the Man- 
ning potential V(0) = -A sech 2 (0/2p) + B sech 4 (0/2p) 
p0[ and the double-Gaussian model V(4>) = l/2(0/p) 2 - 
In cosh(0er//9 2 ), |MJ do not allow exact or quasi-exact 
solvability of the associated Schrodinger equation. We 
note that the double-Morse potential, which arises in sev- 
eral physical contexts such as hydrogen-bonded chains 
p2fl, is closely related to the DSHG potential. 



II. THE DOUBLE SINH-GORDON MODEL 

We start with the continuum representation of a model 
Hamiltonian describing a system capable of undergoing 
a displacive transition: 



H 



dx 
T 



rn 



Ft + V&) + 



(i) 



where I is the lattice spacing, m the mass of particles 
(ions) and Co the velocity of low-amplitude sound waves 
(phonons) in the associated discrete problem [g]. Here 
the potential is 



V DS hg{4>) = (C cosh 20 -n) 2 



(2) 



where £ is a positive parameter. While the value of n is 
not restricted in principle, it has to be a positive inte- 
ger for the QES property to hold, along with n > £ in 
order to have a double well. The two minima occur at 
cosh 20o = n/C with V m m = and a local maximum at 
= with V(0) = (n — C) 2 - We depict the potential for 
a few values of n and £ in Fig. 1. 

The DSHG potential written in the form (||) has all 
the generic features of a double-well potential such as 
Landau-Ginzburg, but allows for much greater analytic 
progress. Below we find exact solutions for (1) a kink, 
(2) a kink lattice, (3) phonon dispersion, and (4) the 
first few eigenvalues and eigenfunctions of the transfer 
operator at certain temperatures both above and below 
the transition, thereby allowing analytic calculation of 
the PDF and correlation functions in the thermodynamic 
limit. 




n=3, zeta=0.05 
n=3, zeta=0.1 



n=2, zeta=0.05 
n=2, z9ta=0.1 



FIG. 1. The DSHG potential for n = 2, 3 and ( = 0.05, 0.1. 



III. KINK AND KINK LATTICE SOLUTIONS 

The single kink solution for the double sinh- Gordon 
potential is nothing but the finite energy solution of 
the equation of motion i<?0 2 (:x;) = Vdshg{4>)- The 
kink/antikink located at xq is given by P,[lO[ 



4>{x) — ± tanh 1 




C , v x - £o 
tann — 
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2(n 2 - C 2 ) 



n > C 



(3) 



Here g 



and tanh </) = \/(n - ()/( n + 0- The 



constant g is often introduced in condensed matter treat- 
ments as a phenomenological parameter and controls 
the kink size. Note that instead of the potential (||) 
one could consider a more general potential of the form 
a(£ cosh 260 — n) 2 . It is easily seen that the kink so- 
lution (j3j) is unchanged except for the replacement of g 
by g/alr. For simplicity, throughout this paper we take 
a, b = 1 since solutions for the general case are trivially 
obtained from here. Traveling kink solutions arc obtained 
by boosting to velocity v via x — > (1 — w 2 ) _1 / 2 (a; — vt). 
The total energy density e(x) = V((f>) + ±0 2 . Thus, the 
energy (or rest mass) of a kink is 



E„ = 



e(x)dx 



= 4^n^n 2 - C 2 tanh -1 



/ n-C 
n + C 



-2£(n 2 -C 2 )- 
(4) 



The topological charge of the kink is 

Q- 



Tr-dx = 0(+oo) — 0(— oo) 
ox 



= 2 tanlT 1 i / - — - = 20o 



(5) 
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The kinetic energy contribution of a moving kink with 
velocity v is obtained from 4*{(x — vt) / ^yl — v 2 ) 



Eki 



m 2 m 2 
ax— d>f = — — v 
2 Yt 2 



2n 



— — ^ tanh 1 



'n-C 



(6) 



To calculate the contribution to the internal energy 
and specific heat from kinks we treat them as a one- 
dimensional gas of indistinguishable, independent, static 
particles of width ~ 2£. The free energy contribution is 
then given by 



F, 



K 



-Nk B T- ( 1 



1 , 27rfc B T 

o ln — 7^2 I ex P 

2 m*B z 



k B TVT 



(7) 



where B is a phase space normalization constant. 

In order to understand kink-antikink interactions, it is 
useful to construct a kink lattice solution which is ob- 
tained by twice integrating the equation of motion 



± 



2 . f*W d<t> 

(x — x ) = / = , 

h{x ) \Jv{4>) - V 



with the boundary conditions for a finite length L — x 2 — 
x%: <t> x {x\) = (f>x{x 2 ) = 0, and V(0(xi)) = V((j){x 2 )) = 
Vq. Here Vq and xq are integration constants and V(x) > 
Vq for x S {x\, x 2 ). The solution is 



{>l(x) = i tanh 1 ^tanh^i 
tanh <j}\ 



x — .To 



k 



tanh ( 



2( sinh (f>i cosh <f> 2 y 2 



(8) 



where d is the periodicity of the kink lattice, K(k) is the 
complete elliptic integral of the first kind with modulus 
k, sn(x, k) is the Jacobi elliptic function, and 



cosh 2^1,2 = ^ T 



k 2 = 



< V < V(cj) = 0) = (C - n) 2 



(9) 



Equivalently, 



tanh 



n+C 



>1,2 = 



(1 + k 2 )C T k' 2 (n + () + \/k ,4 n 2 + 4fc 2 C 2 
(1 + k 2 )C, T fc /2 (n - C) + \/k' 4 n 2 + 4k 2 ( 2 



(10) 



and the characteristic length of a kink in the lattice 



(1 + fc 2 )C 2 - k' 2 {n 2 - C 2 ) + Cy/fc' 4 " 2 + 4fc 2 C 2 
(1 + fc 2 )C 2 + C Vk' 4 n 2 + 4fc 2 C 2 



(11) 



where k! = \/l — k 2 is the complementary modulus. 

The topological charge (per period) in the lattice prob- 
lem Ql = 2(Pl(K) — 2<fii is smaller than the single kink 
case Q = 2(j) = 2tanh~ 1 y/(n — Q)/{n + £). Similarly, 
the kink size in the kink lattice, is also smaller than 
the free kink size £. 

The energy of the kink lattice per period (i.e., the en- 
ergy per kink-antikink pair plus the interaction energy) 
is obtained, after considerable algebra, using 4>l{x) in 
Eq. (1) and evaluating various integrals involving Jacobi 
elliptic functions |23]: 



El = 4£l 



(n + () 2 K + ^(n 2 - ( 2 )(K - E) - 4<n 



where E(k) and n(tanh 2 0i,fc) are complete elliptic in- 
tegrals of second and third kind, respectively. Note that 
in the dilute limit [k — ► 1, d — > oo) the divergences 
in K{k) and n(tanh 2 </>i, k) exactly cancel out and we 
recover the single kink result E s . The interaction en- 
ergy as a function of separation {i.e., k or d) is given by 
Ei n = El — 2E S . 

The repulsive (attractive) kink-kink (kink-antikink) in- 
teraction in the asymptotic limit can be calculated by 
letting k — ► 1 or k' — * and expanding El up to or- 
der fc' 4 . However, for the DSHG model it can directly 
be obtained from the asymptotic form of the single kink 
solution by using Manton's formula [E4[, 



U(r) = ± 



2 n 2 - C 2 



■ exp 



IV. KINK LATTICE THERMODYNAMICS 

One can compute various "thermodynamic" quanti- 
ties associated with the kink lattice similar to the sine- 
Gordon case [p5[. The internal energy per kink U/N = 
El/2 and L/N — d imply that the thermodynamic pres- 
sure 



P 

1 
2 



dU 
~dL 



9 * t 

2n—{n 



T=0 
2 /-2\l/2 



-(n 2 



e)-(n~o 2 

(12) 

where /3 = sin -1 tanh 02- The enthalpy is H = (U + 
PL) = Nfi, with the chemical potential given by 



^ = 2E(n 2 -e)C 



2<z 



-F((3,k)-1 
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where F((3, k) and E(f3, k) are incomplete elliptic inte- 
grals of the first and second kind, respectively. Similarly, 
the isothermal compressibility Ht = — (l/L)(dL/dP)T 
can be calculated: 

Ht = \ {1 + [^/Hd^ L /dk)][E/k' 2 K - 1]} 

x{n/z L )V^em k) - (n 2 - c 2 )(e/e L ) 

-K/fc)^" 2 - ( 2 (df3/dk)(d£, L /dk) 
x{E(/3,k)-F((3,k)}y 1 (13) 

where d(3/dk — cosh 2 (d tanh fo/dk) and d^/dk can 
be obtained from the above expressions for tanh</(>2 and 
£l as a function of k. 

V. STABILITY OF KINKS 

Since the transformation <j) — > tanh connects the 
equation of motion (see below) for the 4 and DSHG 
models, the linear stability of the kink and kink lat- 
tice solutions in the two models is directly related. 
We find that the DSHG kink stability equation is the 
Heun's equation which has at least two bound states: 
(1) the usual zero frequency Goldstone mode and (2) a 
kink shape oscillation mode. Similarly, for linear sta- 
bility of the kink lattice we can write tanh(0£,£) = 
tanh(0i) + c/)l(x, t). Assuming a harmonic time variation 
4>l{x, t) = ip(x) exp(iujt) we find that the kink lattice sta- 
bility equation is also Heun's equation. A special case is 
the Lame equation of order two (y = 2) which arises in 
the context of the stability of the (/> 4 model kink lattice 
solution: 

4>xx + [A(u, k) - v(l + v)k 2 sn 2 x\i) = . 

The solutions of this equation are the Lame functions 
p6[ with 2^+1 eigenvalues. This implies a maximum 
of five bound states. Note that in the integrable case of 
the sine-Gordon equation one gets the Lame equation of 
order one [y = 1) implying a maximum of three bound 
states for the kink lattice and only the Goldstone mode 
for the kink ^5|. However, the DSHG Heun's equation 
does not reduce to a Lame equation. These results apply 
to linear stability (kink-phonon scattering) and cannot 
be readily generalized to the case of nonlinear stability. 

VI. PHONONS 

The equation of motion for the field is given by 

mfitt - 4>" + 2C[Csinh40 - 2nsinh20] = . 

This equation can be linearized around (1) <j) = or (2) 
(j) = 4>o leading to higher and lower energy phonons, re- 
spectively. For phonons around </> = the above equation 



can be approximated by 

16 

rrufa - 4>" + (8C 2 - 8<)0 + — (4( 2 - <)0 3 = . 

For small amplitude oscillations (<fi <C 1) we retain only 
the linear term. Assuming <f){x,t) = 0exp[i(ga; — <-<Jqt)] we 
get the dispersion relation 

lo 2 = q 2 - 8C(n - C) , 

where the frequency uj q will be real only for finite q 2 > 
8C(n— C). Of course, for n < £ (when <fi = is the minium 
of the potential) , the frequency uj q is always real. 

Similarly, the phonon dispersion around the minima 
±<fio for this model is 

u q =q +8{n -C ) ■ 

Note that the frequencies are real for all q > (and 
n > () and the quantity 4(n 2 — C 2 ) is equivalent to the 
coefficient \ A\ of the quadratic term in the corresponding 
(j) 4 model [|| . The associated phonon contribution to the 
free energy is 

with S being the lattice constant and (3 = l/fc^T. The 
contribution to the internal energy and specific heat from 
phonons and kinks can be directly calculated using F V ib 
and Fk, respectively. However, here we do not consider 
the contribution to the free energy due to kink-phonon 
interactions p9[ . 

VII. TRANSFER OPERATOR FORMALISM 

The transfer operator method transforms the problem 
of finding the canonical partition function for a system to 
an exactly equivalent problem of finding the eigenvalues 
of a certain integral operator. In the continuum limit, 
this problem can be further reduced to finding the en- 
ergy eigenvalues of a related quantum mechanical prob- 
lem (i.e., a Schrodinger equation with an effective poten- 
tial). In addition, the transfer operator formalism also 
provides exact expressions for correlation functions. The 
key property of the DSHG potential ^ is that the n low- 
est eigenvalues and wavefunctions are known exactly, in 
contrast to other QES systems. This enables exact eval- 
uation of the PDF and the correlation functions C\ and 
Ci (see below) at several temperatures. Aside from val- 
idation of Langevin results against these quantities, the 
exact single kink solution can be checked against that 
obtained from the low temperature Langevin field con- 
figuration. 

Turning now to the computation of i? c z, we note that 
this calculation can be divided into two parts: a triv- 
ial Gaussian integration over the field momentum, and 
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a computation of the configurational partition function, 
which via the transfer integral method becomes equiv- 
alent to solving a Schrodinger-like equation The 
Hamiltonian for the DSHG theory is 



H = dx 



— IT 

2 



Vdshg(4>) 



(15) 



and this leads to the Schrodinger-like equation for the 
eigenvalues and eigenfunctions of the transfer operator, 



'WW 



■* fc + (C cosh 20- n)H k 



(16) 



This is a QES system for which, at 2(3 2 = 1, using results 
for a related potential from Ref. ||] , the eigenstates of the 
first n levels can be found for n = 1, 2, 3, 4. (We have been 
able to extend this to the cases n = 5, 6.) However, what 
one really wants is to consider a given fixed-n theory and 
obtain eigenstates at different temperatures. It is easy to 
see from Eq. (|l^) , by simple rescaling, that solutions of a 
fixed-n theory at certain values of /3 are the same as the 
solutions of another theory (different n and () at 2/3 2 = 1. 
Depending on the chosen value of n, exact solutions are 
available at different fixed values of (3. Here, we restrict 
ourselves to one such family (n — 2) which allows the 
exact computation of the first few eigenstates at 8/3 2 = 
m 2 (m = l,---,6). For illustration, three examples of 
the (unnormalized) ground states are given below (see 
also Fig. 2). The first (high temperature, /3 2 = 1/8) has 
an eigenfunction with a single peak while the other two 
(lower temperatures, 1 = 1/2 and 9/8) have a double 
peak: 



*O(0)| /3 2 = 1 



exp ( — -£ cosh 20 



cosh exp ( — cosh 20 



3C+ (1 + V1 + 9C 2 ) cosh 20 



x exp ( cosh 20 



(17) 



with corresponding ground state energies, Eq — 4 + £ 2 , 
E = C 2 - 2( + 3 and E = ( 2 + f - § y/T+W respec- 
tively. The PDF for the classical field is just the square 
of the normalized ground state eigenfunctions. Solutions 
at higher energies and other values of f} are given in Ap- 
pendix A. 

Once the eigenvalues of the transfer operator are 
known, they can be used to compute the correlation func- 
tions d{x) = (0(O)0(x)) andC 2 = (<50 2 (O)<50 2 (a;)), using 

C 1 (x)=^K* fe |0|*o)| 2 exp[-/3|a ; |(^-£;o)] , (18) 

k 

C 2 {x) =X)K*fc N 2 |*o) I 2 exp [-(3\x\(E k - E )] . (19) 



It is apparent that at large distances, C\ and C2 are 
dominated by the lowest state with nonvanishing matrix 
elements: the first excited state in the case of C\ and the 
second excited state in the case of Ci . Since Eg , E\ , and 
E2 are known at certain temperatures, the large distance 
behavior of these correlation functions can be found ex- 
actly and compared with the results from simulations. 
Static structure factors may also be calculated in much 
the same way. 




FIG. 2. PDFs at three values of [3 2 with the continuum 
exact solutions shown as solid lines. Results from Langevin 
simulations are superimposed at the /3 values: 1/8 (triangles), 
1/2 (diamonds), and 9/8 (squares). 

At this point, it is important to mention the connection 
between the "quantum" calculations and kink physics. 
In the context of kink statistical mechanics, it is usual 
to introduce a phenomenological description of kinks as 
particles in a grand canonical ensemble. However, this is 
unnecessary, and all such thermodynamical information 
can be extracted directly from the Schrodinger descrip- 
tion of the transfer operator. For example, the kink den- 
sity has been obtained in this way in Ref. p7| . Simpler 
quantities like C\ and C2 have obvious natural interpre- 
tations in terms of kinks. The C\ correlation length is 
related to the kink/antikink spacing and increases mono- 
tonically as (3 increases (Fig. 3). The behavior of the C2 
correlation length is more subtle, since C2 is not sen- 
sitive to domain size. At both high temperatures (no 
kinks) and low temperatures (number of kinks exponen- 
tially suppressed), the correlation length is dominated by 
the thermal phonon contribution. However, in the range 
of temperatures close to the kink transition the fluctua- 
tions on the kink length scale become important and can 
dominate Ci . At these temperatures one might expect a 
maximum in the C2 correlation length, on the order of 
the kink size, and this is indeed what we observe numeri- 
cally (Fig. 4). The Schottky anomaly in the specific heat 
PjJ arises for the very same reason. 
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FIG. 3. Comparison of the inverse correlation lengths from 
Ci(x) obtained via Langevin simulations [/3 2 : 1/2 (diamonds), 
and 9/8 (squares)] and from the large \x\ continuum exact 
results (solid lines). 
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FIG. 4. The numerical and exact inverse correlation 
lengths from C2(x) are compared for the same three tem- 
peratures and with the same conventions as in Fig. 2. The 
solid line is the large \x\ continuum exact result for f3 2 =9/8. 
The largest correlation length is at the intermediate value of 
P (see text). 

VIII. QES SOLUTIONS 

The explicit QES solutions (exact eigenvalues and 
eigenf unctions) for 2/3 2 = 1 and n = 1,2,3,4 are given 
by Razavy || while we have been able to extend these 
results to n — 5, 6. We now explain in more detail how 
to relate these results to a theory with a given value of 
n (say n = 2) but with different values of f3. We start 
from the Schrodinger equation (|l^) at 2/3 2 = 1 and arbi- 
trary n (= 1, 2, • • • , 6). On multiplying both sides of this 
equation by 4/n 2 , it immediately follows that the exact 
solutions for a given value of n, £, E at 2/3 2 = 1 are equiv- 
alent to the solution of the DSHG equation at n = 2 but 



at 8(3 2 — n 2 , E — AE/n 2 , ( = %/n. Clearly, instead of 
fixing n — 2, one could choose any of the six values of n 
and for each case one knows the exact solutions for that 
theory at six different temperatures. As an illustration, 
the explicit QES solutions (exact eigenvalues and eigen- 
functions) for n = 1, 2 at six different temperatures are 
presented in Appendix A. 

IX. CONNECTION WITH 4 AND DOUBLE 
SINE-GORDON MODELS 

As a final point, we consider the relationship of the 
DSHG theory to the more familiar Landau-Ginzburg 
model. Scrutiny of the static equation of motion [<fi 2 = 
V((f))] reveals an important connection between the kink 
(and kink lattice) solutions of the </> 4 model and the dou- 
ble sine-Gordon (DSG) and DSHG models. Consider the 

</> 4 potential V A (u) = [(n + Qu 2 - (n - 0] 2 . The substi- 
tution u — tanh cj> takes the (static) equations of motion 
over to the DSHG equations. The alternative substitu- 
tion u = tan</> leads to a DSG model. This means that 
all known solutions of the </> 4 equation of motion can be 
directly taken over to the DSHG and DSG equations of 
motion (and vice versa). As one use of this interesting 
relationship, the DSG kink lattice solution (not known 
heretofore in the literature) can be written down directly 
in case Vdsg = (C c °s20 — n) 2 : 

4> L =±tan~ 1 ^tan^i sn ^~ X ° ,fc^ , (20) 

simply by using the substitution tanh (f> — > tan <p in the 
kink lattice solution of the DSHG model. 

This connection enables us to write down by inspection 
not just the kink solutions but their total energy as well, 
which is often a very tedious task. Similarly, knowing 
the linear stability of the 4> A kink and kink lattice, im- 
portant results follow for the stability of corresponding 
solutions of the DSG and DSHG models (see Sec. V). 
Moreover, since we know that the DSHG model is an ex- 
ample of a QES system, and considering the very similar 
way in which the DSG and DSHG models are related to 
</> 4 , it is logical to conjecture that some DSG model may 
also be a QES system. Indeed, this is the case, and we 
have found several exact eigenvalues and cigcnfunctions 
at many temperatures for a particular DSG system. The 
exact statistical mechanical results for that model, simi- 
lar to the DSHG results presented here, will be reported 
elsewhere 

Alternatively, if we start from the </> 4 potential 14 (u) = 

\(n + C) ?i2 + ( n ~~ C)l 2 ; the transformation u = tan0 
takes the (static) equations of motion over to the model 
Vdsg — [C cos 20 — n] 2 . For n < ( we get the doubly 
periodic DSG potential with two types of kinks ( "small" 
and "large"). When n > £, we get the singly periodic 
DSG potential with the 27r-kink solution (see Fig. 1 of 
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Ref. p^l). The associated Va{u) in this case is a sin- 
gle well potential. However, the real transformation does 
not lead to a QES DSG system unlike the transformation 
above. 

We note that the transformations u = tan <p and 
u = tanh cj) connect the equation of motion of the </> 
model to that of the DSG and DSHG systems, respec- 
tively. On the other hand, the transformations u = cos (j> 
and u = cosh <f> connect the equation of motion of the (f> 4 
model to that of the exactly solvable sine-Gordon (Vsg = 
sin 2 (ft) and sinh-Gordon (Vshg — sinh 2 </>) systems, re- 
spectively, provided we start with Va(u) — (it 2 — l) 2 . We 
also note that the number of periodicity in the trigono- 
metric potentials (DG, DSG, triple sine-Gordon (TSG), 
• • •) corresponds to the number of wells in the associated 
hyperbolic potentials (SHG, DSHG, triple sinh-Gordon 
(TSHG), ■ • •). We conjecture that the equations of mo- 
tion for the TSG and TSHG models can be directly con- 
nected with that of the (f> 6 model via appropriate trans- 
formations. In addition, while the TSG equation [ p^|j29| 
is known to arise in nonlinear optics the TSHG poten- 
tial [e.g. Vtshg = A cosh 2cj) + B cosh 4(f) + C cosh 6(f)] 
can serve as a model system to study first-order phase 
transitions. 



X. SIMULATION RESULTS 

The numerical study of kink statistical mechanics was 
carried out via Langevin simulations. The advantage of 
such simulations is that representative field configura- 
tions (see Fig. 5) are available and can be analyzed to 
compute such quantities as the kink density, finite tem- 
perature kink profile, the spatial kink distribution, kink 
transport and nucleation, etc. This sort of analysis is not 
possible with Monte Carlo techniques. 




100000 120000 140000 160000 180000 200000 



FIG. 5. Sample kink configurations from a lOOif section 
of a 400K unit simulation with f3 2 = 9/8. 

The additive noise Langevin equation for the double 
sinh-Gordon model is 

d% (f> = dl x (j> - r)d t (f> - 4((C cosh 2(f) - n) sinh 2(f) + F(x, t) 



where the viscosity r\ and the Gaussian white noise F are 
related by the fluctuation-dissipation theorem: 

^F(x,t)F(x',t')^ = 2 V (3- I 6{x - x')5{t - t') . 

The lattice versions of the above continuous equations 
(with periodic boundary conditions) were solved using 
standard techniques [j30| . Random initial conditions were 
driven to equilibrium and the results sampled in time 
thereafter to yield time averaged PDFs, etc. The use 
of the Langevin technique for obtaining thermodynamic 
quantities is straightforward and as shown below, can be 
remarkably accurate. The availability of exact solutions 
is very useful since they set stringent criteria for accuracy 
that must be met by numerical methods. In the simu- 
lations reported below we used a lattice size of typically 
512K points with a lattice spacing A = 0.025. The time 
step was taken to be e = 0.005. [We employed both ^Je 
(Euler) and e 2 (Runge-Kutta) order stochastic integra- 
tion algorithms.] These values are substantially smaller 
than what has been the norm so far in Langevin simula- 
tions of other 1 + 1 dimensional field theories |HJ. In fact, 
based on our experience with the double sinh-Gordon sys- 
tem it is likely that errors in previous Langevin analyses 
have been as high as 30%. 

Fig. 2 shows the striking agreement between the nu- 
merically obtained and the exact continuum PDFs at 
three temperatures: The worst case departure is at the 
level of parts per thousand. The comparisons for the 
inverse correlation lengths are given in Figs. 3 and 4. 
DSHG system parameters are n = 2, £ = 0.05. For 
Ci, the numerical values are 1/A = 0.1425 (f3 2 = 1/2) 
and 1/A = 0.012 (f3 2 = 9/8) as compared to the exact 
values in the continuum theory of 0.14142 and .0105, re- 
spectively. The small offset between the continuum and 
lattice calculations is due to the finite value of the lattice 
constant and is consistent with estimates from higher- 
order contributions to the transfer integral J^[. The kink 
number density at low temperatures is related to the cor- 
relation length via Nk — 1/4A J2^J . This general relation- 
ship is borne out in the DSHG simulations. 

The high quality of these numerical simulations im- 
plies that the PDF can now be used directly to compute 
thermodynamic quantities at any temperature. Since 
the PDF is just the square of the ground state wave 
function of the Schrodinger equation, one can use it to 
compute the ground state energy Eq numerically, from 
which the internal energy (U = dEo/d/3), the free energy 
(F = E /f3), and the entropy (S = f3dE /df3 - E ) can 
all be computed in a straightforward way. The specific 
heat involves two (3 derivatives and is difficult to obtain 
with good accuracy but in this case, the standard energy 
fluctuation method is quite effective. The use of the PDF 
complements traditional techniques utilizing energy fluc- 
tuations in Langevin simulations which are not suited to 
free energy and entropy calculations. 

The QES nature of the DSHG theory allows not only 
the exact computation of E at several temperatures, but 
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also of dEo/d(3, using first order perturbation theory: 
dE o /dj3\ 0=0o = (^o,d 2 ^ /d(j) 2 ) where f3 is one of the 
special checkpoint temperatures. Thus the internal en- 
ergy U and the entropy S can also be found exactly at 
these temperatures. Once again, these quantities can be 
used to validate numerical work over a broad range of 
temperatures. 



XI. CONCLUSION 

In this paper we have analytically obtained the free en- 
ergy of the kink-bearing DSHG theory in 1 + 1 dimensions 
at several temperatures. To our knowledge, this is the 
first such analytical calculation at several temperatures. 
We have compared the analytical results against large 
scale Langevin simulations and found excellent agree- 
ment between the two. In view of this agreement, the 
study of PDF-based thermodynamics via Langevin sim- 
ulations appears quite tractable. It is worth emphasizing 
that the calibration of Langevin simulations constitutes 
a fundamental problem in the numerical study of nonlin- 
ear dynamical systems. Numerical convergence of results 
is certainly an important test that does not require the 
knowledge of an exact solution. However, standard fi- 
nite difference error analysis cannot be applied directly to 
Langevin simulations because of their stochastic nature. 
Moreover, the error and convergence analysis for nonlin- 
ear stochastic PDEs with spatio-temporal noise remains 
to be completely worked out. For this reason, a nontriv- 
ial test system for which analytic solutions are available 
at several temperatures provides a valuable tool for test- 
ing the Langevin code. In fact based on our experience 
in the DSHG case, it is likely that previous Langevin- 
based analyses have errors as high as 30%. We believe 
that apart from the intrinsic physical interest of the the 
double-well DSHG model, this system has the potential 
of becoming a standard benchmark problem in numerical 
simulations of Langevin systems. 

We have also analytically studied the static equations 
of motion of the DSHG system and found an interesting 
connection between the 4 -theory, DSHG theory and a 
DSG theory. As a result, the kink and lattice kink solu- 
tions as well as kink total energy can be directly written 
down for any one of them in terms of the other. This is 
interesting as the lattice kink solution of this DSG case, 
which was unknown so far, can be immediately written 
down. Further, because of the intimate connection be- 
tween the DSHG and the DSG case, it is reasonable to 
expect that a DSG theory may also be a QES system. 
We have been able to obtain several eigenstates of the 
DSG case at various temperatures from this identifica- 
tion j3^| . It is clearly of great interest to carry out high 
quality Langevin simulations of the DSG case, which is 
a periodic system, and compare the agreement between 
the analytical and the simulation results. We hope to 
report on this problem in the near future. 
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APPENDIX A: QES SOLUTIONS 

We obtain two sets of exact solutions for the 
Schrodinger equation (|l6|). 

Set-I (n = 1): Note that eo, as it occurs in the free 
energy, is related to the ground state energy Eq by 
e a = E - V min . For n = 1, V min — for ( < 1 while 
V m i n — 2£ — £ 2 — 1 for £ > 1. However, in either case, 

E n E m . 



(1) For 2(3 2 = 1: 
tf o (0) = exp -| cosh20 , E = l + C 2 

(2) For 2/3 2 = 4 



(Al) 



#oW>) = cosh 0exp[-£ cosh 2cj>] , E = C? - ( ■ 
#1(0) = sinh exp[-C cosh 20] , E x = £ 2 + C - 

Note that ei - e = E t - E = 2(. 
(3) For 2/3 2 = 9: 



4 ' 
3 

4 ' 

(A2) 



4r o (0) = [6C + (1 + V1 + 36C 2 ) cosh 20] 



x exp 



3C , „ / 

cosh 20 

2 



Eo 



C 2 -^VTT36?+9, 



\&i(0) = sinh 20 exp 
5 



-— cosh 20 
2 Y 



Ei=( z 



9 ' 



* 2 (0) = [6C- (V1 + 36C 2 - l)cosh20] 
3C 



x exp 



■ cosh 20 



E 2 = C 2 + i(2Vl + 36C 2 + 7) 



(A3) 
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Note that ei - e = flV 1 + 36<C 2 - 1) and e 2 - e 
IV1 + 36C 2 - 

(4) For 2/3 2 = 16: 



*o(0) = [12Ccosh?!.+ (2-4C + 2V1 -4C+ 16C 2 ) 
x cosh 30] exp[— 2£ cosh 20] , 

E ° = 1^( 16 C 2 ~ Vl - 4C + 16C 2 - 8( + 11) , 
lb 

* i (</>) = [12C sinh + (2 + 4C + 2^1 + 4( + 16C 2 ) 
x sinh 30] exp[— 2£ cosh 20] , 
1 



E * = ttt(16C 2 + 8C - Ay/l + 4C + 16C 2 + 11) , 
16 

* 2 (0) = [12Ccosh0 + (2 - 4C - 2y/l - 4C + 16C 2 ) 
x cosh 30] exp[— 2£cosh20] , 

£ 2 = ^(16C 2 - 8C + Ay/l - 4C + 16C 2 + 11) . 
16 



(A4) 



Note that ei - e = C + |\/l-4C+16C 2 - 
iV! + 4 C+16C 2 and e 2 - e = ±0 -4C + 16C 2 . 

Set-II (n = 2): In this case Kim for £ > 2 is 4£ - ( 2 - 4 
while V^ mm = in case £ < 2. 



(1) For 2/3 2 = 1: 
*o(0) = cosh0exp 
5'i(0) = sinh0exp 



- — cosh 20 



-— cosh 20 



E = C - 2C + 3 , 

Si = C 2 + 2C + 3 . 

(A5) 



Note that t\ - e Q = 4£. 
(2) For 2/3 2 = 4: 



4- (0) = [3C cosh + ( 1 - C + V1-2C + 4C 2 ) cosh 3C] 
x exp[— C cosh 20] , 

Eo = |(4C 2 - Vl - 2C + 4C 2 - 4C + 11) , 

4-1(0) = [3C sinh + (1 + C + V1 + 2C + 4C 2 ) sinh 30] 
x exp[— C cosh 20] , 

Ei = J(4C 2 + 4C - Ay/l + 2( + 4C 2 + 11) , 



4f 2 (0) = [3Ccosh0+ (1 - C - y/l - 2C + 4C 2 )cosh30] 
x exp[— C cosh 20] , 
1 



E 2 = -(4C 2 - 4C + 4^1 - 2C + 4C 2 + 11) . (A6) 



Note that e 1 - e = 2C, + y/l-2Q + 4( 2 - y/l + 2Q + 4<C 2 
while e 2 - e = 2^1 - 2C + 4( 2 . 



(3) For 2/3 2 = 1/4: 



*o(0) = exp 



C 



cosh 20 



E = C + 4 



(A7) 



(4) For 2/3 2 = 9/4: 



*o(0) = [3C + (1 + V1 + 9C 2 ) cosh 20] 
3C 



x exp 



■ cosh 20 



^o = C 2 -^v^W+f 



^i(0) = sinh 20 exp 
Et-C +-, 



— ^ cosh 20 



* 2 (0) = [3C - (V1 + 9C 2 - 1) cosh20] 

K 

X exp — — cosh 20 , 
( 2 + ^(2y/TT9e + 7) . 



Eo 



(A8) 



[y/l + 9C 2 - 1] and e 2 



Note that e% 

We now show that the first five and six levels of the 
DSHG potential as given by Eq. (||) can be written down 
exactly in case n — 5 and n — 6 respectively (at 2(3 2 = 
1). Using the trivial scaling, one can then obtain exact 
solutions for say n = 2 but at 2j3 2 = 25/4 and 20 1 = 9. 
For example, for n = 5, Eq^a are solutions of the cubic 
equation 

E(E + A)(E + 16) - 64C 2 £- 768C 2 = , (A9) 
and the corresponding eigenfunctions are given by 



0(0) = v(4>) exp 



— — cosh 20 



(A10) 



where 



^0,2,4 = 6C — E cosh 20 



2(E 



E+16 



cosh 40 



E 



1,3 



-10t2 v /4C 2 + 9 , 



, 3 = 4C sinh 20 + (3 ± y/AC, 2 + 9) sinh 40 



(All) 



Note that here E = E - C, 2 - 25. 

For n=6, Bq 2 4 are solutions of the cubic equation 



(J5 + 9)(E + 25)(E + 6C + 1) - 52( 2 E 
-820C 2 - 120C 3 = , 



(A12) 



and 
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wo, 2, 4 = 8Ccosh</> - (E + 6(+ 1) cosh 30 
+ 2 C (£ + 6C + l) cosh 
£ + 25 

Note that £1,3,5 are solutions of the same cubic equation 
with £ replaced by — £, and 

Ulj3i5 = 8Csinh0- (£-6C + l)sinh30 
+ 2 C (£-6 C+ l) sinh 
£ + 25 

Here £ = £ - C 2 - 36. 



[1] D. J. Scalapino, M. Sears, and R. A. Ferrell, Phys. Rev. B 
6, 3409 (1972). 

[2] J. A. Krumhansl and J. R. Schrieffer, Phys. Rev. B 11, 
3535 (1975). 

[3] A. G. Ushveridze, Quasi- exactly Solvable Models in Quan- 
tum Mechanics (IOPP, Bristol 1994); F. Cooper, A. Khare, 
and U. P. Sukhatme, Phys. Rep. 251, 267 (1995). 

[4] A. Saxena and S. Habib, Physica D 107, 338 (1997). 

[5] A. Khare, S. Habib, and A. Saxena, Phys. Rev. Lett. 79, 
3797 (1997). 

[6] A. D. Bruce, Adv. Phys. 29, 111 (1980). 

[7] A. D. Bruce, T. Schneider, and E. Stoll, Phys. Rev. Lett. 
43, 1284 (1979); A. D. Bruce, in Nonlinear Phenomena 
at Phase Transitions and Instabilities, edited by T. Riste 
(Plenum, New York, 1981); K. Binder, Z. Phys. B 43, 119 
(1981). 

[8] S. Habib and G. D. Lythe, (in preparation). 
[9] M. Razavy, Phys. Lett. A72, 89 (1979); Am. J. Phys. 48, 
285 (1980). 

[10] S. N. Behera and A. Khare, J. Physique C6, 314 (1981). 
[11] K. B. Joseph and B. V. Baby, J. Phys. A 16, 2685 (1983); 

B. Dey, J. Phys. A 19, 1735 (1986). 
[12] P. T. Dinda, Phys. Rev. B 46, 12012 (1992). 
[13] A. M. Dikande and T. C. Kofane, Solid State Commun. 

89, 283 (1994); ibid 559 (1994). 
[14] H. Konwent, Phys. Lett. A 118, 467 (1986); H. Konwent, 

P. Machnikowski, and A. Radosz, J. Phys. A 28, 3757 

(1995). 

[15] R. K. Bullough, P. J. Caudrey, and H. M. Gibbs in Solitons, 
eds. R. K. Bullough and P. J. Caudrey (Springer- Verlag, 
Berlin, 1980), Chapt. 3; M. J. Ablowitz, M. D. Kruskal, 
and J. F. Ladik, SIAM J. Appl. Math. 36, 428 (1979). 

[16] E. L. Ince, Proc. R. Soc. (Edinburgh) 45, 102 (1925). 

[17] M. J. Ablowitz, D. J. Kaup, A. C. Newell, and H. Segur, 
Stud. Appl. Math. 53, 249 (1974); H. P. Mckean, Comm. 
Pure Appl. Math. 34, 197 (1981). 

[18] S. E. Trullinger and R. M. DeLeonardis, Phys. Rev. B 20, 
2225 (1979). 

[19] J. F. Currie J. A. Krumhansl, A. R. Bishop, and S. E. 

Trullinger, Phys. Rev. B 22, 477 (1980). 
[20] M. F. Manning, J. Chem. Phys. 3, 136 (1935). 



[21] G. A. Baker, Jr. and A. R. Bishop, J. Phys. A 15, L201 

(1982) ; G. A. Baker, Jr., A. R. Bishop, K. Fesser, P. D. 
Beale, and J. A. Krumhansl, Phys. Rev. B 26, 2596 (1982); 
J. H. Chen, M. E. Fisher, and B. G. Nickel, Phys. Rev. 
Lett. 48, 630 (1982); P. D. Beale, Phys. Rev. B 27, 5804 

(1983) . 

[22] E. Matsushita and T. Matsubara, Prog. Theor. Phys. 67, 

1 (1982); ibid 71, 235 (1984); E. S. Kryachko, Chem. Phys. 

143, 359 (1990); O. Yanovitskii, G. Vlastou-Tsinganos, 

and N. Flytzanis, Phys. Rev. B 48, 12645 (1993). 
[23] P. F. Byrd and M. D. Friedman, Handbook of Elliptic 

Integrals for Engineers and Scientists, 2nd ed. (Springer, 

Berlin, 1971). 
[24] N. S. Manton, Nucl. Phys. B 150, 397 (1979). 
[25] N. Gupta and B. Sutherland, Phys. Rev. A 14, 1790 (1976). 
[26] B. Horovitz, Phys. Rev. B 35, 734 (1987); A. Saxena and 

W. Cao, Phys. Rev. B 38, 7664 (1988). 
[27] F. J. Alexander, S. Habib, and A. Kovner, Phys. Rev. E 

48, 4284 (1993). 
[28] D. K. Campbell, M. Peyrard, and P. Sodano, Physica D 

19, 165 (1986); 

[29] R. M. Gundersen, Intern. J. Nonlin. Mech. 26, 389 (1991). 
[30] A. Griner, W. Strittmatter, and J. Honerkamp, J. Stat. 

Phys. 51, 95 (1988). 
[31] F. J. Alexander and S. Habib, Phys. Rev. Lett. 71, 955 

(1993), and references therein. 
[32] S. Habib, A. Khare, and A. Saxena (in preparation). 



10 



